#Dependencies
library(tidyverse)
library(haven)
library(readxl)
library(knitr)
library(janitor)
library(stringr)

#Retrieve the categories of different scenarios from metadata
codes <- read_excel("~/Dropbox/AR6/AR6_metadata.xlsx",
                    sheet="meta_Ch3vetted_withclimate") %>%
    select(Scenario,Category) %>% 
    distinct(Scenario,.keep_all=TRUE)

ar6 <- read.csv("~/Dropbox/AR6/AR6_scenarios.csv") %>% 
    left_join(codes) %>% #Add the categories to our main data
    filter(Category %in% c("C1","C2","C3"))

rm(codes)


#BECCS ANALYSIS

beccs <- ar6 %>% 
  filter(Variable=="Carbon Sequestration|CCS|Biomass") %>%
  arrange(X2100) %>%
  filter(!is.na(X2100) & !is.na(X2050)) %>%
  group_by(Category) %>%
  select(-c(Variable,Unit,Region))

options(digits=2)
tab_BECCS <- beccs %>% summarise(LOW_exceeds_low_2050 = round(mean(X2050 > 700)*100),
                            LOW_exceeds_med_2050 = round(mean(X2050 > 1400)*100),
                            LOW_exceeds_high_2050 = round(mean(X2050 > 3300)*100),
                            LOW_exceeds_low_2100 = round(mean(X2100 > 700)*100),
                            LOW_exceeds_med_2100 = round(mean(X2100 > 1400)*100),
                            LOW_exceeds_high_2100 = round(mean(X2100 > 3300)*100),
                            MED_exceeds_low_2050 = round(mean(X2050 > 1200)*100),
                            MED_exceeds_med_2050 = round(mean(X2050 > 2900)*100),
                            MED_exceeds_high_2050 = round(mean(X2050 > 7800)*100),
                            MED_exceeds_low_2100 = round(mean(X2100 > 1200)*100),
                            MED_exceeds_med_2100 = round(mean(X2100 > 2900)*100),
                            MED_exceeds_high_2100 = round(mean(X2100 > 7800)*100))

#AFFORESTATION/REFORESTATION ANALYSIS

aff <- ar6 %>%
  filter(Variable %in% c("Carbon Sequestration|Land Use")) %>%
  arrange(X2100) %>%
  filter(!is.na(X2100) & !is.na(X2050)) %>%
  group_by(Category) %>% 
  select(-c(Variable,Unit,Region))


tab_AR <- aff %>% summarise(exceeds_low_2050 = round(mean(X2050 > 1300)*100),
                            exceeds_med_2050 = round(mean(X2050 > 3800)*100),
                            exceeds_high_2050 = round(mean(X2050 > 6300)*100),
                            exceeds_low_2100 = round(mean(X2100 > 1300)*100),
                            exceeds_med_2100 = round(mean(X2100 > 3800)*100),
                            exceeds_high_2100 = round(mean(X2100 > 6300)*100))